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Cosmic-ray muons and especially their secondaries break apart nuclei (“spallation”) and pro¬ 
duce fast neutrons and beta-decay isotopes, which are backgrounds for low-energy experiments. In 
Super-Kamiokande, these beta decays are the dominant background in 6-18 MeV, relevant for solar 
neutrinos and the diffuse supernova neutrino background. In a previous paper, we showed that 
these spallation isotopes are produced primarily in showers, instead of in isolation. This explains an 
empirical spatial correlation between a peak in the muon Cherenkov light profile and the spallation 
decay, which Super-Kamiokande used to develop a new spallation cut. However, the muon light pro¬ 
files that Super-Kamiokande measured are grossly inconsistent with shower physics. We show how 
to resolve this discrepancy and how to reconstruct accurate profiles of muons and their showers from 
their Cherenkov light. We propose a new spallation cut based on these improved profiles and quan¬ 
tify its effects. Our results can significantly benefit low-energy studies in Super-Kamiokande, and 
will be especially important for detectors at shallower depths, like the proposed Hyper-Kamiokande. 


I. INTRODUCTION 

Astrophysical neutrinos can reveal the extreme physi¬ 
cal conditions in their sources as well as new information 
about neutrino properties. In the MeV energy range, the 
key targets are solar neutrinos and the diffuse supernova 
neutrino background (DSNB). Solar neutrinos have been 
detected for half a century, yet there are still unanswered 
questions m- The upper limit on the DSNB flux is 
within a factor of a few of theoretical predictions Effig. 

Super-Kamiokande (Super-K) is a 50-kton water 
Cherenkov neutrino detector [14], ITS] . Due to its large 
volume, low backgrounds, and long running time, Super- 
K has the best sensitivity to the high-energy, low-flux 
branches of the solar neutrinos and to the DSNB. 

In Super-K, these measurements are background lim¬ 
ited. The dominant background in 6-18 MeV is the spal¬ 
lation background Emun, which consists of beta de¬ 
cays from unstable isotopes produced by muons and espe¬ 
cially their secondary particles M- These backgrounds 
are reduced by associating them with their muon par¬ 
ents, which can be difficult because some of these iso¬ 
topes have long lifetimes (several seconds) compared to 
the muon rate (~ 2 Hz) [T7|. 

In our first paper in this series P9. the only theoreti¬ 
cal study of spallation in water, we calculated the average 
spallation yields in Super-K. We compared the aggregate 
time profile and energy spectrum of spallation decays to 
Super-K measurements, finding agreement within uncer¬ 
tainties. We showed that almost all isotopes are made 
by secondary particles, e.g., neutrons, pions, and gamma 
rays, instead of primary muons. 

In our second paper [50], we showed for the first time 
that almost all spallation isotopes are made in muon- 
induced showers. These showers have high densities 
of secondary particles; they extend only ~ 5 m along 
muon tracks, while the height of the Super-K detector 


is ~ 40 nr. Because showers can be detected through 
their Cherenkov light, this provides a new way to iden¬ 
tify where a spallation isotope might be produced along 
the muon track. 

Earlier, Super-K empirically found variations in the 
Cherenkov light intensity along muon tracks, and a corre¬ 
lation between the position of the peak and the spallation 
decay [5]. They developed a new spallation cut based on 
the measured correlation. Using this, they lowered the 
analysis energy threshold for the DSNB search. However, 
the physical cause of the light variations and their corre¬ 
lation with spallation were unexplained. Also, they did 
not apply this cut to their solar neutrino analysis. 

Our finding that most spallation isotopes are made 
in showers explains Super-K observations, except one. 
Their reconstructed muon light profiles are much broader 
and have much smaller amplitude than those expected 
from showers. Here we show how this discrepancy can 
be explained by shortcomings of the Super-K reconstruc¬ 
tion method, and how to improve it. We explore applica¬ 
tions of better-reconstructed profiles. Our results should 
greatly benefit their solar neutrino and DSNB analyses. 
Although we use Super-K as an example and attempt to 
model its main present features, our focus is more gen¬ 
eral. 

For our calculations, we use the simulation package 
FLUKA (version 2011.2c.0) 51, 55], It incorporates all 
the relevant physics for muon interactions in water. Our 
physics choices for FLUKA are the same as in our pre¬ 
vious papers p2i [20]. We simulate throughgoing muons 
vertically down the center of the Super-K detector; our 
results can be applied to more general cases. The muon 
spectrum is shown in Fig. 1 of Ref. [T9] . At the Super-K 
depth (2700 meter water equivalent), the average muon 
energy is 270 GeV pH? ] [55 ] . 

One difference in our setup here is that our simulation 
region is the whole Super-K inner detector (ID), whereas 
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in our previous papers we used only the fiducial volume. 
The ID is a cylinder 33.8 m in diameter and 36.2 m in 
height [13] , It is separated from the outer detector by 
opaque walls (including ceiling and floor), where photo¬ 
multiplier tubes (PMTs) are mounted [13]. The fiducial 
volume is an analysis region inside and smaller than the 
ID (22.5 kton versus 32 kton) 13] . The PMTs collect 
light emitted in the whole ID, so we use it for our simu¬ 
lation volume. 

This paper is organized as follows. In Sec. [jl] we dis¬ 
cuss the basics of shower physics and muon light profiles. 
In Sec. m we review the Super-K reconstruction method 
and how to improve it. In Sec. IV we explore further 
applications of better-reconstructed shower profiles and 
quantify how much they could improve the spallation cut. 
We conclude in Sec. El 


II. MUON CHERENKOV LIGHT PROFILES 

Relativistic charged particles in water emit Cherenkov 
light along their paths. The Cherenkov photons propa¬ 
gate through the detector, occasionally getting scattered 
or absorbed. Some of the photons reach PMTs and are 
detected. The light intensity (number of photons per dis¬ 
tance) emitted by a singly charged particle per unit dis¬ 
tance is constant, independent of the particle type and 
energy [23] . The total number of photons is proportional 
to the energy deposited, and their arrival positions and 
times carry information about the event geometry. 

When cosmic-ray muons pass through Super-K, they 
produce charged secondary particles, such as electrons 
and pions. (With these generic terms, we typically mean 
e ± and 7r , we separate 7r°.) The production and energy 
loss of secondary particles is prompt, much faster than 
muons crossing the Super-K ID (~ 100 ns). It is thus 
not straightforward to separate the Cherenkov light from 
cosmic-ray muons and their secondaries. 

We call the Cherenkov light intensity along a muon 
track the muon light profile. Its fluctuations reveal 
secondary production, because the light intensity from 
muons is constant. To better describe the production of 
secondaries, we separate it into two steps: a primary 
muon directly produces daughter particles, and these 
daughter particles subsequently produce other secondary 
particles. 

We quantify charged particles not by their number, but 
by the distance they travel, which is proportional to their 
Cherenkov light emission. In FLUKA, charged particles 
are propagated by track segments. The Cherenkov light 
emission from each segment is proportional to its length. 


A. Muon-produced Charged Particles 


The energy loss rate of a muon is 



a + /3E^ t , 


( 1 ) 


where the brackets indicate averaging over distance 1251 . 
The a term is for ionization loss, and the fdE^ term is for 
radiative loss. At the Super-K depth, where ( E = 270 
GeV, a ~ 2.9 MeV cm -1 and P{E^) ~ 0.7 MeV cm -1 [23, 
[ 251 . 

Ionization is muons losing energy by scattering bound 
electrons [25]. We can further divide this term based on 
the energy of the outgoing electrons. When their energy 
is small, this is restricted ionization; when it is large, 
this is delta-ray production |25j . Restricted ionization 
loss is a continuous process, while delta-ray production 
is discrete interactions. The boundary between these two 
cases is somewhat arbitrary [22]; we set it to be the elec¬ 
tron Cherenkov threshold (kinetic energy 0.257 MeV). 
The average total ionization energy loss for a muon that 
travels vertically through the ID is ~ 10 GeV, with ~ 6 
GeV due to the restricted loss and ~ 4 GeV due to delta- 
ray production [19] . 

Radiative processes include pair-production, 
bremsstrahlung, and photonuclear interactions 125] . 
These are muons interacting with nuclei and producing 
electron-positron pairs, gamma rays, pions, and other 
mesons. All of these processes have a large energy 
transfer for each discrete interaction (up to hundreds of 
GeV [25]), and the interaction rates are low. The total 
radiative energy loss through the ID is ~ 3 GeV. All 
the energy that goes into radiative processes is carried 
by secondary particles, and mostly dissipates through 
ionization of these secondary particles. Because of this, 
there is a near-constant relationship between energy loss 
and Cherenkov yield. For the production spectra of 
daughter particles by muons in Super-K, see Fig. 7 of 
Ref. [20]. 

The fluctuation levels of the energy losses determine 
the features of muon light profiles. The restricted ioniza¬ 
tion loss has negligible fluctuations along muon tracks, 
and among different muons. The delta-ray production 
and radiative energy losses have large fluctuations. Even 
though, on average, these terms are smaller than the re¬ 
stricted ionization, they can be much larger for individ¬ 
ual muons. We show how these features affect muon light 
profiles in the next subsection. 

The energetic daughter particles from delta-ray pro¬ 
duction and radiative processes produce many secondary 
particles by inducing electromagnetic and hadronic show¬ 
ers [2D] . A shower is a series of repetitive interactions 
where particles interact and multiply in number and de¬ 
crease in energy. When the average particle energy in 
the shower is too low to create new particles, particles 
range out by ionization. An electromagnetic shower is 
mostly gamma rays producing electrons and positrons 
by pair production, and electrons and positrons produc¬ 
ing gamma rays by bremsstrahlung. A hadronic shower 
is mostly charged pions producing multiple charged and 
neutral pions, and neutral pions decaying to gamma rays 
and inducing electromagnetic showers. 

An important energy scale for showers is the critical 
energy E c \ showers develop when the average particle 


dx 




3 




FIG. 1. An example of a real (simulated) muon light profile. 
The height is scaled to the light from a single muon. The 
thin red line at 16 m indicates the peak position of the light 
profile. 


FIG. 2. An example of a reconstructed muon light profile from 
Super-K (data from Fig. 2 of Ref. [5]). The height is scaled 
to the light from a single muon. The thin red line indicates 
the peak position of the light profile. Compared to the profile 
in Fig. [l] this peak is much wider and much shorter. 


energy is above it, and die out below it. In water, E c ~ 
100 MeV for electromagnetic showers |25j and E c ~ 1 
GeV for hadronic showers [20]. Most showers in Super-K 
have energies ~ 1-300 GeV [20] , 

Because most shower particles are energetic and have 
small deflections, showers look like long thin cylinders in 
real space (an example is shown in Fig. 1 of Ref. HO])- In 
terms of Cherenkov light production, the dominant con¬ 
tribution for either kind of shower comes from electrons 
near E c ~ 100 MeV. The average deflection of electrons 
is (cos 9 Z ) ~ 0.8. Showers in Super-K extend ~ 5 m in 
the longitudinal direction, and ~ 10 cm in the lateral 
direction (in hadronic showers, this can be ~ 1 m) m- 
For a detailed discussion of shower geometry in Super-K, 
see Sec. IIC of Ref. [20], 

Muon daughter particles with energy below E c do not 
induce showers. The most common such particles are 
electrons. Even though these electrons do not induce 
showers, they are important for our discussions because 
they emit Cherenkov light. We refer to them as low- 
energy delta rays, i.e., electrons with energy above their 
Cherenkov threshold but below about 100 MeV. We in¬ 
clude in this a ~ 5% contribution of low-energy electrons 
plus positrons from pair production. 

Shower physics is the key to understanding how to use 
Cherenkov light to tag spallation backgrounds. There 
are abundant electrons in electromagnetic and hadronic 
showers, so showers can be observed through their 
Cherenkov light. There are few (many) pions and neu¬ 
trons in electromagnetic (hadronic) showers, and they 


efficiently make isotopes. By observing the light from a 
shower, Super-K could identify its position and thus the 
position where spallation decays might occur [2D]. How¬ 
ever, the light profiles observed by Super-K [3123 look 
very different from what we expect from showers. We 
take a closer look at this next. 


B. Real vs. Reconstructed Muon Light Profiles 

To utilize the correlation between showers and spalla¬ 
tion decays, we need to know what showers look like in 
the detector. In this subsection, we first study real muon 
light profiles, focusing on the shower shape. By real we 
mean what reconstructed profiles would look like if every 
detected photon were reconstructed to its correct emis¬ 
sion position; the profiles are simulated. We then look 
at the differences between these profiles and those recon¬ 
structed from Super-K data. 

Figure [l] shows an example of a real muon light profile. 
We simulate a vertical throughgoing muon in the ID and 
plot the total path length of relativistic charged particles 
relative to the muon path length. This is equivalent to 
the Cherenkov light intensity in units of that from a sin¬ 
gle muon. We adopt this unit because it is not affected 
by experimental effects, such as photon absorption and 
detector efficiency. For this and similar figures, we use a 
bin size of 0.5 m, to be consistent with Super-K |9]. 

The area under the curve is proportional to the muon 
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energy loss. A height of 1 corresponds to only restricted 
muon ionization energy loss. Any height larger than 1 is 
due to ionization of additional charged particles, which 
are produced through delta-ray production or radiative 
processes. 

Though just a single example, Fig. [T] is representative. 
When muons are not showering, their profiles all look 
similar, with the relative intensity fluctuating between f 
and 2 due to the muon plus low-energy delta rays (the 
height is sometimes < 1 due to binning issues). The 
average level is about 1.3, which corresponds to an extra 
energy loss of 0.5 MeV cm -1 due to low-energy delta 
rays. The muon profile in Fig. [l] shows one energetic 
shower of energy ~ 15 GeV. (This example peaks near 
the center of the detector; showers can occur anywhere 
along the muon track.) It is quite typical, extending ~ 5 
m along the muon track and with a height of about 30. 
For a shower of energy Eq, the peak height is typically 
2-3 (So/GeV) [23. 

The distance between the peak position of the light 
profile and a candidate signal event is used to determine 
the probability of the event being a spallation decay [9]. 
In other words, Super-K keeps the peak position and dis¬ 
cards information about the shape of the muon light pro¬ 
file. In Sec. |IV[ we discuss how this can be improved. 

The number of showers and the shower energies vary a 
lot for each muon [2D] ■ The more energetic the daughter 
particle, the more rare it is. Low-energy showers are thus 
more common than high-energy showers. The average 
number of showers above 1 GeV per vertically through- 
going muon in Super-K is 0.4. It is most common to have 
zero or one shower per muon. Multiple showers along one 
muon track happen less than 10% of the time. 

Super-K measured the variations of muon light inten¬ 
sity in Ref. [D[. Though their approach is fairly general, 
it is not based on showers, in which electron deflections 
play a crucial role. In the next section, we discuss how 
this affects muon light profile reconstruction. 

Figure [2] shows an example of a reconstructed muon 
profile from Super-K (another is Fig. 4.2 of Ref. f27j). 
Their original y axis is in units of number of photoelec¬ 
trons detected. We convert that to our units by using 
7p.e. ~ 1 MeV [25], which includes Cherenkov photon 
yield, photon absorption, and detector efficiency, etc. 
This relation is not exact, but it does not affect our dis¬ 
cussions. 

This light profile varies along the entire muon track, 
showing a prominent peak in the middle. The full width 
of the peak is about 20 m and its height is about 10. We 
estimate that the excess light corresponds to about 15 
GeV, similar to the example shown in Fig. [l] Beyond 32 
m, the falloff in intensity is probably because this muon 
track left the ID. 

It is puzzling why the reconstructed profile looks like 
this. Even though the shower in Fig. jT] is only one re¬ 
alization, it is representative of a shower of a similar 
energy. The success of the current Super-K cut indi¬ 
cates that their average reconstructed peak position, at 


least for large muon energy losses, is quite robust. In the 
next section, we explain the differences in the muon light 
profiles, and how to improve the Super-K reconstruction 
method to get accurate light profiles. 

III. SHOWER RECONSTRUCTION 

It is important to understand the discrepancy between 
the real and reconstructed muon light profiles. First, it 
is crucial to the proof that showers are the true cause 
of the variation in muon light intensity and that spalla¬ 
tion isotopes are made in showers. Second, as we pointed 
out in Ref. [20], the correlation function Super-K mea¬ 
sured is not as sharp as the one we calculated using real 
shower profiles. If Super-K could better reconstruct ev¬ 
ery shower, even if its energy is low or if it is accompanied 
by other showers, it would improve the efficiency of the 
spallation cut. 

It is difficult to reconstruct the muon light profile. Dur¬ 
ing a shower, there are many charged particles emitting 
light at nearly the same time and position, but pointing 
in different directions, so it would be very difficult to re¬ 
solve individual rings in the pattern of PMT hits on the 
wall. Furthermore, the light from showers at specific lo¬ 
cations must be separated from the continuous light from 
the muon itself and from low-energy delta rays. 

In this section, we use our knowledge of showers, which 
are not mentioned in the Super-K paper, to examine 
the Super-K reconstruction method and its results. We 
start by reviewing the setup of the equation that Super- 
K used for reconstruction. Then we study in detail its 
properties and how its solutions are affected by physical 
and detector limitations. We present ways to improve 
their method, demonstrating that we can reconstruct the 
muon light profile with high fidelity, including identifying 
showers and measuring their energy, position, and extent. 

We first assume that the times of individual detected 
photons at a given PMT can be measured separately. 
Then, in Sec. |III D[ we discuss complications to that in 
Super-K and possible solutions. 

A. Super-K Reconstruction Method 

To reconstruct a muon light profile, Super-K performed 
backward fitting using individual PMT hits, without con¬ 
sideration of their correlations. For each PMT hit (taken 
here to be one detected photon), they measured its posi¬ 
tion and time and solved for where along the muon track 
the photon was emitted. They repeated this for all hits 
from one muon, and got the number of detected photons 
as a function of the emission position along the muon 
track, i.e., the muon light profile. 

There is no confusion about which muon to associate 
a photon with, because photons from one muon arrive at 
the walls within ~ 100 ns, and the average time between 
muons is ~ 0.5 s (muon bundles are treated separately). 
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FIG. 3. Diagram for reconstruction of the muon Cherenkov 
light profile. A photon is emitted from the muon track jl at 
time f e and distance xi, propagates a distance X 2 , and hits 
a PMT at time t and position r. The blue triangle (with 
corners marked by to,t e ,t) is described by Eq. § and the 
green curve by Eq. j3j; they cross at the solutions of Eq. 0. 
of which only one is marked. 


In addition, the short duration of the signal means that 
photons from unrelated low-energy decay backgrounds 
can be ignored. Finally, because the water is so clear, the 
effects of light absorption and scattering are minimal. 

Figure [3] illustrates the geometry. A muon enters the 
ID at time to and moves along the direction (1 with speed 
c. A Cherenkov photon is emitted at time t e , a distance 
x\ = c(t e — to) along the muon track by either the pri¬ 
mary muon or a secondary charged particle. The photon 
propagates a distance x 2 to the PMT with group veloc¬ 
ity c/rig [25H5T] (n g = 1.38 for water). The photon hits 
the PMT at time t, where its position is r relative to the 
muon entry point. Because this method treats one de¬ 
tected photon at a time, in essence redrawing Fig. [3] for 
each one, it is easy to accommodate muons at arbitrary 
positions and angles, as well as PMT hits on the ceiling 
or floor of the ID. The azimuthal angle of each PMT hit 
in fixed detector coordinate is needed to define the plane 
of Fig. [3j after which it is not used. 

In Fig. [3j the muon position at every instant, and thus 
ft and to, is known because the entry and exit points 
and times are determined using inner and outer detector 
information. The PMT hit time t and its position relative 
to the beginning of the muon track, r, are measured. The 
angle 9 is known immediately. 

The angle a can be obtained after solving for x± and 
t e . For photons emitted by primary muons, a is the 
Cherenkov angle ao, defined by the photon phase velocity 


via cosao = 1 /n p h [251431] . Its value in water is ao = 
42°, with n p h = 1.33. Notice here that n g and n p h have 
similar values but are not equal. 

The distance X\ and time t e can be calculated as the 
joint solutions of two separate constraints, as illustrated 
in Fig. [3] The spatial constraint is 

x\ = r 2 + x\ — 2a;ircos 9 , (2) 

which is satisfied by all points along the solid blue line of 
the muon track. The time constraint is 


x\ + n g x 2 = c(t — t 0 ), 


(3) 


which is satisfied by all points along the green curve for 
which the total time — accumulated as a muon for a 
distance X\ and as a photon for a distance x 2 — is t — to- 
If n g = 1, this curve would be an ellipse with foci at the 
points associated with t 0 and t ; instead, it is a fourth- 
order polynomial, and we show a relevant section. 

The joint solutions are defined by a quadratic equation 
in x±, obtained by combining the constraints: 


4- 1 


xf+2 


r-fi- 


z{t - to) 


Xi- 


\t-t 0 f 


= 0 . 


(4) 

This is the same as Eq. (4.3) in Ref. [27] : we provide 
more details on the origin and solutions of this equation. 
In principle, there could be zero, one, or two solutions. 
Because every observed photon is emitted at some point 
on or close to the muon track, there should always be 
at least one real solution. In Sec. EH we show that 
numerical issues can make it appear that there are none. 
The Super-K procedure is to measure the coefficients in 
this equation and solve for X\, keeping all solutions that 
are real and fall into a reasonable range. 

The Super-K reconstruction method is compatible 
with the properties of showers, because it is approxi¬ 
mately true that all the light is emitted from a point 
moving along a straight line. The longitudinal spread 
of particles at one instant of a shower is ~ 0.1 m. For 
an electromagnetic shower, the transverse extent is also 
~ 0.1 m; for a hadronic shower, it can be ~ 1 m. Finally, 
typical muon track deflections are small, coincidentally 
~ 0.1 m over the height of the ID. As we discuss be¬ 
low, the typical resolution for shower reconstruction is of 
order a few meters, so these effects are negligible. 

This method is powerful because it accounts for all of 
the observed Cherenkov light. For relativistic muons, a 
photon is emitted at a fixed angle a = a 0 relative to the 
muon track; in this case, the emission point of the pho¬ 
ton could also be obtained from a simpler linear equation, 
e.g., Eq. ([5]). (An earlier Super-K paper [32j, on atmo¬ 
spheric neutrinos, reconstructed the muon light profile 
assuming that all emission was at the Cherenkov angle 
relative to the muon.) However, for electrons, which can 
be significantly deflected, the angle a varies, and the full 
quadratic equation is needed. The major shortcoming of 
this method, in either the quadratic or linear case, is that 
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it neglects correlations between PMT hits, as discussed 
in detail in Sec. mu 

As we have shown, the Super-K reconstructed light 
profiles are inconsistent with what we expect from a 
shower (or any known process). Yet, at first glance, the 
Super-K reconstruction method looks correct. To under¬ 
stand the differences between the real and reconstructed 
light profiles, we take a closer look at the nature of solu¬ 
tions of Eq. Q. 


B. Understanding the Reconstruction Solutions 

To demonstrate how reconstruction works, we first sim¬ 
ulate the Cherenkov light pattern that Super-K observes 
on the ID walls from muons and their secondaries. We 
define the z axis to point in the direction along the muon 
track. For every PMT hit on the walls, Super-K can mea¬ 
sure its position z = r- (i and time t. By light pattern, we 
mean the ( z , t) plane filled by all the PMT hits from one 
muon and its secondaries. The azimuthal angle for each 
PMT hit defines a plane for Fig. [3j when the light pattern 
is constructed, those azimuthal angles are discarded. 

To calculate the light pattern, we make some rea¬ 
sonable simplifications; a full study by Super-K will be 
needed to fine-tune the details. We use a cylinder with 
the same radius (16.9 m) as the Super-K ID, but take 
it to be much taller, so that all light is collected on the 
side walls, instead of the ceiling or floor. We take the 
muons to be vertically downgoing at the center of the 
cylinder, so that there is azimuthal symmetry on aver¬ 
age. For more general cases where the muon is tilted or 
shifted from the center, the appearance of Fig. [4] changes 
and we discuss it separately. We use (m, ns) as the units, 
suppressing their display below. The muon enters the ID 
at (0, 0). 

We generate and propagate Cherenkov photons geo¬ 
metrically with our own codes, ignoring light absorp¬ 
tion and scattering. For each track segment (~ 1 cm) 
in FLUKA charged particle propagation, we effectively 
propagate light with intensity equal to its length along 
its Cherenkov cone from the midpoint of the segment. 

The Super-K data are discrete in 2 (all points on the 
surface of a PMT are taken to be the center of the PMT) 
but continuous in t (though there is smearing due to time 
resolution). For computational reasons, it is simpler for 
us to take z to be near-continuous (bins of width 0.05 m) 
and t to be discrete (bins of duration 3 ns). For each z 
value, we calculate t and round it to the nearest bin. If 
a segment is tilted away from the vertical, we uniformly 
distribute its light between the minimum and maximum 
z. The discreteness of 3 ns in t, which is comparable to 
the timing resolution of the Super-K PMTs, is roughly 
equivalent to 1 m in z. We discuss resolution further 
below. 

As a check, we generate the light patterns with the 
Cherenkov light propagation in FLUKA, including light 
absorption and scattering, and the results are consistent. 



Arrival position [m] 


10 1 

10 ° 

10 1 

io -2 

10 -3 

io ~ 4 


FIG. 4. The time and position (z = r-pt) pattern of Cherenkov 
light on the ID walls for one example of a muon and its sec¬ 
ondaries. We simulate an infinitely tall cylinder, such that all 
the photons hit the side walls instead of the ceiling or floor. 
The arrival position range is thus larger than the height of 
the Super-K ID (36.2 m). The intensity scale is approximate; 
the light from muons is about 1 on the scale. The real muon 
light profile for this example is shown in Figs. [5] and [6] 


Furthermore, when we use the Super-K reconstruction 
method, we recover light profiles similar to theirs. 

Figure [4] is an example of the light pattern produced 
by a muon and its secondaries. The diagonal band with 
the highest intensity is due to the muon and the most 
forward electrons in one moderate shower, a couple of 
small showers, and some low-energy delta rays (which 
have (cos 6 Z ) ~ 0.85). The remainder of the intensity, 
less than 20%, which arrives at given positions at later 
times, is due to significantly deflected electrons; that in 
the bottom left corner of the figure is due to very de¬ 
flected electrons that produce light near the ceiling of 
the detector, where the muon cannot. 

How much of this Cherenkov intensity (per bin of time 
and position) can be detected? The first impression 
shows blocks of size 1 m, though a zoom-in reveals the 
0.05-m bins; both have height 3 ns. Because 1 m of a 
charged particle corresponds to ~ 10 3 PMT hits, Super- 
K should be able to detect some light from such blocks 
with intensities > 10 -3 (or from 0.05-m bins with inten¬ 
sity > 10 —2 ), assuming each is integrated over 3 ns. 

The pattern in Fig. El from a muon (or any track seg¬ 
ment along the 2 axis) can be understood from simple 
arguments. The setup is shown in Fig. [3j but here we 
assume a vertical central muon. The photons emitted at 
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a point Xi have 


z = x i + 16.9 cot ao = x\ + 19, 


and 


t = 



16.9 rig 
sin «o 


/c 


ad 

c 


118, 


(5) 

( 6 ) 


where 16.9 (m) is the Super-K ID radius and ao is the 
Cherenkov angle. (For noncentral or nonvertical muons, 
these and the following expressions can be generalized 
by changing 16.9 to rsind.) The offsets correspond to 
the position shift and time delay for light to reach the 
walls after the muon first enters the ID. These equations 
combine to give the pattern 


z = ct — 16.9(7-^-cotao) (7) 

\ sin a 0 ) 

in Fig. 0 This line is broadened to the left, into a band, 
due to finite detector time and position resolution. 

The pattern in Fig. [4] from deflected electrons also has 
a characteristic shape. For deflection by an angle 9 Z , the 
minimum and maximum z for the propagated photons 
are X\ + 16.9cot(ao ± 0 Z ). In a shower, there are many 
electrons in a short distance, some with very large deflec¬ 
tions. For light emitted at all angles from a single point 
x\, the complete pattern is a hyperbola, 

c ! (*-t £ f = _ 16 92i 

U 9 

placed to the left of the line defined by the muon. The 
lowest point of the hyperbola, zi OW) comes from light that 
travels perpendicular to the wall, which reveals the point 
of emission through X\ = Zi ow . The right-hand side of the 
hyperbola is populated by light from electrons aligned 
close to the muon track; the left-hand side by light that 
is moving upward in the detector, due to very deflected 
electrons. Because shower electrons are forward-peaked 
((cos$z) ~ 0.8), the intensity on the hyperbola falls as 
the electron deflection increases. Outside the range of 
Fig-0 the Cherenkov intensity is nonzero but negligible. 

For other muon positions or orientations, the appear¬ 
ance of Fig. 0 changes. For a noncentral but vertical 
muon, the line from the muon light turns into a band. In 
this case, photons emitted at the same point but at dif¬ 
ferent azimuthal angles travel different distances to the 
wall. Consequently, for a fixed z value on the wall, it 
gets hit by photons emitted at different positions along 
the muon track, and they take different times to reach z. 
For nonvertical muons, the width of the band varies be¬ 
cause a nonvertical muon can be considered to be many 
small tracks of noncentral vertical muons. These effects 
do not change our results because our reconstructions are 
based on Fig. 0and Eq. ( 0 , which are fully general for 
each PMT hit. We obtain good results for other muon 
positions and orientations. In effect, for each range of 


azimuthal angle and height in detector coordinates, the 
pattern of PMT hits looks like that shown in Fig. 0 
The key to understanding Fig. 0 and hence the re¬ 
construction method — is electron deflections. We quan¬ 
tify these by the discriminant of Eq. 0 , 


4 x 16.9 2 
sin 2 a 



(9) 


which determines the nature of the solutions (for central 
vertical muons). If the phase and group velocities were 
identical, A would be zero for photons emitted by muons. 
For these, A ~ 1, much smaller than typical values for 
electrons. To simplify the discussion, we approximate 
A ~ 0 for muon light. 

Most electrons are quite forward (a ~ ao), so A is 
usually small; because measured PMT hits correspond 
to physical solutions, A must be positive. When an elec¬ 
tron is aligned with the muon, A = 0 and the two solu¬ 
tions merge, corresponding in Fig. 0 to the green curve 
from Eq. ( 0 ) being tangential to the blue muon line from 
Eq. @. When an electron is deflected, A > 0 and there 
are two real solutions, as in the example shown in Fig. 0 
As A increases, the two real solutions split further apart. 
Although both are physical, one is correct and one is not, 
and these cannot be distinguished on an event-by-event 
basis. The Super-K reconstruction method keeps both 
solutions if they are within a reasonable range, which 
leads to a problem of overcounting PMT hits (about 35% 
in the example of Fig. 0. 

The typically small discriminant amplifies the effects 
of the detector time and position resolution. The mea¬ 
sured A can become negative, corresponding to no real 
solutions, due to shifts in the measured quantities used 
in the coefficients of Eq. 0 . In Fig. 0 the green curve 
would be slightly displaced from the blue muon line, with 
no crossings; in Fig. 0 the PMT hits would be slightly to 
the right of the line defined by the muon. The Super-K 
reconstruction method discards such cases, which leads 
to a problem of undercounting PMT hits (about 40% in 
the example of Fig. 0 . 

These numerical problems also mean that true solu¬ 
tions depend sensitively on the measured values of (z, t). 
In Fig. 0 the near-straightness of the green curve means 
that slight movements of it or the blue muon line lead to 
large changes in the solutions. In Fig. 0 it is difficult to 
separate hyperbolas with different emission points X\ by 
looking at their right-hand sides, where the Cherenkov 
intensity is greatest. More quantitatively, 


xi(z,t) 


xi(zo,t Q ) + 


2.6 

71 


(ct — z)(Sz 


cSt), 


( 10 ) 


where Sz and St describe how incorrect the values of 
(z,t) are. When A is small, the error term scales as 
~ 20/\/A ~ 0.5| cos a—cos ao| _1 m, which can be several 
meters for typical shower electrons. This is significantly 
larger than the PMT position or time resolution because 


of the near-cancellation in A. (As we show in Sec. IIIC 








when A is near zero, we can solve a linear equation in¬ 
stead of the quadratic.) Importantly, this tells us that 
the best reconstructions come from the worst electrons , 
i.e., those with the largest deflections and A values. 

The nature and precision of the solutions can also be 
affected by the presence of hadronic showers. For these, 
transverse displacements of some shower particles can be 
~ 1 m from the muon track, especially for large show¬ 
ers. This means that A values are shifted compared to 
the case with no deflection. This can be seen from the 
fact that there are more negative A values and they can 
have larger absolute value, compared to photon hits from 
electromagnetic showers. The effect is small enough that 
we can ignore it here, but large enough that it could help 
identify hadronic showers, which are rare but which pro¬ 
duce nearly all isotopes. 

With a better understanding of the solutions, we now 
have insights as to why the Super-K reconstructed profile 
looks very different from real shower profiles. We next 
consider how to improve their method. 



C. Improving the Reconstruction Method 


The exploration in Sec. Ill B reveals how the Super-K 


reconstruction method works, as well as three improv¬ 
able limitations. First, when A is large, taking both so¬ 
lutions includes wrong information and overcounts PMT 
hits. Second, when A is negative due to detector reso¬ 
lution, taking zero solutions ignores correct information 
and undercounts PMT hits. Third, when A is small, the 
sensitivity of the solutions to detector resolution dilutes 
better information. 

These limitations result in reconstructed muon light 
profiles with distorted shapes and inaccurate shower en¬ 
ergies. When the shower energy is small, the current 
method might not be able to localize the shower. Multi¬ 
ple showers cannot be resolved either. 

Our goal for improving the reconstruction method is to 
get an accurate muon light profile. This includes locat¬ 
ing the correct shower peak position, getting the correct 
shower shape and shower energy, and resolving multiple 
showers. We improve the Super-K method by addressing 
its three limitations. First, when A is large, we show how 
to select the better solution (Improvement 1). Second, 
when A is negative, we show how to repair it and recover 
a solution (Improvement 2). Third, when A is small, we 
show that, though these solutions help reconstruct the 
complete light profile, it is best to set them aside when 
defining the showers (Improvement 3). 

Figure [5] shows (in gray shade) the Cherenkov light 
profile from a simulated muon and its secondaries. Using 
this specific example, we calculated the (z, t ) data shown 
in Fig. [4j here we use that data as if they were observed, 
attempting to reconstruct an accurate light profile from 
it. 

In this example, the total muon energy loss is 11 GeV. 
There is a medium-sized shower (about 4 GeV, as can be 


FIG. 5. Shower reconstruction with refinement of real solu¬ 
tions. The gray shaded shape is the real (simulated) light 
profile used to produce the example light pattern in Fig. [4] 
The black line is the result of the reconstruction using the 
Super-K method. For the red line, we keep only the better 
one of the two solutions. 


seen from its area) located near 10 m, a smaller shower 
near 5 m, and possibly some smaller ones further along 
the muon track. These are all quite typical in appear¬ 
ance, with the smaller showers being harder to recognize. 
These particular showers are electromagnetic; there are 
harmless isotopes produced by gamma rays at around 5 
m and 30 nr. As explained before, the light outside the 
shower regions is from the muon and low-energy delta 
rays. It may look like there are larger fluctuations in the 
muon and delta-ray light than in Fig. [T] but this is only 
because the overall y scale is smaller due to this shower 
being smaller. 

We choose this example because the biggest shower has 
only moderate energy and because there are two showers 
close to each other. It is a good test of how well the 
reconstruction works, both in terms of getting the correct 
shape of the largest shower and of resolving the small 
showers. 

Figure [5] also shows (black line) the result obtained 
when we use the Super-K reconstruction method. The 
largest shower is found at the right position, which is 
why the new Super-K spallation cut works, as shown in 
Fig. 12 in Ref. (2Dj. However, the shape of this shower is 
badly smeared. The area under the black line is compa¬ 
rable to that in the gray shade, but this is an accident, 
because the Super-K method overcounts and undercounts 
by roughly equal amounts, as noted. The small showers 
are not resolved. 

Figure [5] also shows (red line) the result if we improve 














9 


the Super-K method by keeping only one solution when 
there are two (Improvement 1). We first run the Super- 
K reconstruction method, and record the peak position 
from the reconstruction. We then run the reconstruc¬ 
tion a second time, keeping only the solutions closer to 
the peak. The red line agrees with the black line near the 
peak, as expected, but is lower elsewhere because it is not 
including wrong solutions and overcounting. Though it 
defines the showers better, it is still significantly broader 
than it should be. We note that when the Super-K recon¬ 
struction method produces a wrong peak, due to multiple 
showers or for small showers, Improvement 1 can rein¬ 
force it, but does not cause the problem. Improvement 
3 can fix the problem because “fake” showers would not 
have many deflected electrons, as shown in the Appendix. 

Next we consider the data for which A is slightly nega¬ 
tive due to detector resolution (Improvement 2). Typical 
values are at worst —30 m 2 (nominal resolution ~ 3 m 
[Eq. ©]), corresponding to the level expected from de¬ 
tector resolution. We reconstruct these cases by setting 
A = 0, as would be appropriate for light from the muon 
or very forward electrons. Because the angle of the par¬ 
ticle is known, the quadratic equation reduces to a linear 
equation, e.g., Eq. ([5|. We use the linear equation only 
when A = 0 (or is reset to be), as the quadratic equa¬ 
tion is less sensitive to numerical problems from detector 
resolution except for the smallest positive values of A. 

Figure [6] shows (green line) the result when we also 
recover these formerly discarded solutions (all of our im¬ 
provements of the reconstructed muon profile are cumu¬ 
lative.) This profile is a better match in the peak and 
especially the baseline to the input in the gray shade 
than even the red line in Fig. [5j There is no longer un¬ 
dercounting of light from undeflected particles, such as 
the muon itself. Importantly, the area under the green 
line now matches that in the gray shade. However, the 
shower peak is still too wide, and no secondary showers 
are identifiable. Because the green line traces the non¬ 
showering part of the light profile well, it can be used to 
estimate the energy of the largest shower from its area 
above the baseline; at this stage of refinement, imperfect 
precision increases the width of the shower and decreases 
its height, but conserves its area. 

Finally, we focus on the photons that provide the most 
precise information on the positions of showers (Improve¬ 
ment 3). From Eq. (10), these are the ones with large A. 
We choose A > 100 nr, corresponding to a nominal reso¬ 
lution of 2 m [Eq. ©]■ This comes at a price of keeping 
only ~ 10% of the PMT hits, corresponding to ~ 1 GeV. 
We solve the quadratic equation, keeping only the solu¬ 
tion closest to the shower peak as determined with the 
Super-K method. We correct the normalization of the 
blue line by adding a constant baseline of 1 for a muon 
and by setting the shower energy above the baseline to 
match that of the green line. 

Figure [6] also shows (blue line) the result obtained us¬ 
ing only the most deflected electrons. The agreement 
with the input shown in the gray shade is excellent. Com- 



FIG. 6. Shower reconstruction with two more improvements. 
The gray shaded shape is the real (simulated) light profile. For 
the green line, we repair unphysical data to recover muonlike 
solutions that were previously discarded. For the blue line, we 
select only the light from electrons with large deflections to 
focus on reconstructing only the showers; we add 1 to account 
for the muon light. 


pared to previous results, it is much narrower, localizing 
showers better. Only this method clearly defines multi¬ 
ple showers. We added in the muon baseline to facilitate 
comparison, but the underlying method ignores the light 
from the muon and most of its low-energy delta rays. 
That is, it focuses on the light in showers, where nearly 
all the isotopes are made. 

We show all three improvements in this order to best 
explain the physics. In practice, the first step is Im¬ 
provement 3, which is to pick out the most deflected elec¬ 
trons. Next is Improvement 1, which is to select only one 
solution for each PMT hit (for the most deflected elec¬ 
trons). Lastly, one can follow Improvement 2 to get the 
correct total energy, then scale the profile from Improve¬ 
ment 3+1 to the correct energy. However, one can also 
skip Improvement 2, and get the total energy simply by 
counting the total number of PMT hits, then scaling the 
profile to the correct energy. 

We have demonstrated that the Super-K reconstruc¬ 
tion technique can be significantly improved, leading to 
an accurate muon light profile, even when there are mul¬ 
tiple showers. However, there are some complications re¬ 
garding practical implementation, which we discuss next. 
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D. Towards Practical Implementation 

So far we have assumed that Super-K can measure the 
position and time of each detected photon within the 
precision noted in Sec. EH However, this is not always 
possible with the present electronics. Here we discuss 
the implications and how to achieve the aim of shower 
reconstruction anyway. 

A PMT hit is the basic observable and must be defined 
carefully. In Super-K, it is defined by the total number of 
detected photons within a time window and the time of 
just the first photon. The number of detected photons is 
determined by the accumulated charge of the photoelec¬ 
trons produced. In Super-K, the time window is ~ 400 ns 
(Michael Wilking, private communication; Michael Smy, 
private communication). As assumed above, the position 
and number of detected photons are well defined, but the 
individual times are not, which reduces the available in¬ 
formation. 

This effect is important for reconstructing muon- 
induced showers in Super-K. The light yield of a vertical 
throughgoing muon is high, corresponding to several de¬ 
tected photons per PMT, and more if there are large 
showers. As shown in Fig. [4j the light from the muon al¬ 
ways arrives at a given PMT before that from a shower. 
Because of this, most PMTs lose the timing informa¬ 
tion on the light from showers. Despite this, Super-K 
found reasonable reconstructions in Ref. 9]. where they 
weighted the solutions corresponding to each PMT by 
the total number of detected photons. 

Much of the data needed for reconstruction are not 
affected by this limitation. The key is to identify cases 
where light from the muon does not reach the PMTs. 
The most significant reason is due to geometry. For ver¬ 
tically downgoing muons, the most common case, their 
Cherenkov light cannot reach the PMTs in roughly the 
top half of the detector. This can be seen in Fig. [4j 
the height of the ID is 36.2 m, and the muon light begins 
only at a depth of 19 m. We emphasize that the PMTs in 
the top half of Super-K can detect photons from show¬ 
ers anywhere in the detector. Indeed, the further the 
direction of the shower light is away from that of the 
muon light, the better. Another reason is due to fluctua¬ 
tions. Some PMTs that could have been reached by the 
muon Cherenkov light will not be triggered, and these 
will properly register late-arriving light. 

To check the effects of the timing limitations, we con¬ 
structed a second simulation, which is a more faithful 
representation of Super-K. (We do not use this simula¬ 
tion for our main results because it complicates the dis¬ 
cussion of the underlying physics.) The simulated region 
follows the true Super-K ID geometry with the ceiling 
and the floor. Individual PMTs are mounted on the ID 
walls with realistic sizes and spacing. We record each 
photon hit with its total charge and first-hit time, as op¬ 
posed to treating photoelectrons individually. For the 
reconstruction, we repeat the Super-K method and our 
improvements. For Improvement 3, not only do we se¬ 


lect the photons with large deflection, but also hits on 
the ceiling and in the top half of the detector, where no 
muon light is expected. Our reconstructed profiles rea¬ 
sonably trace the true profiles and can pick out showers 
occurring even near the bottom of the detector. Thus we 
are confident that the properties of the PMT electronics 
will not significantly affect our results. 

Longer term, the ability to reconstruct showers could 
be improved by installing new electronics that allow for 
pulse-shape discrimination, or at least enough informa¬ 
tion to separate detected photons that arrive a few tens of 
ns apart. New technologies for Cherenkov light detection 
with excellent position and time resolution are extremely 
promising for improving shower reconstruction [338435] . 

In the near term, the most promising possibility is to 
go beyond the framework of the Super-K reconstruction 
method and our improvements, and take advantage of 
the ideas proposed in the following subsection. 

E. Towards Better Reconstruction Methods 

The Super-K reconstruction method, including our im¬ 
provements, works reasonably well. However, it has fun¬ 
damental shortcomings. It neglects the correlations be¬ 
tween different photons from the same charged particle, 
i.e., the Cherenkov ring pattern. It neglects the correla¬ 
tions between different electrons emitted from the same 
position, i.e., the shower angular distribution. And it 
neglects the correlations between electrons emitted from 
different positions, i.e., the shower longitudinal profile. 

Better methods should be possible. Here we sketch 
three promising ideas, each for a different energy range; it 
may be possible to combine them. A good reconstruction 
needs only to provide the number of relativistic charged 
particles accompanying the muon, and some information 
about their angular distribution, each in bins of size ~ 0.5 
m along the muon track. The muon and the shower each 
produce a lot of light, ~ 7000 PMT hits per GeV, which 
provides a lot of information for such modest goals. We 
have had encouraging conversations with Super-K col¬ 
laborators about specific codes that could be adapted to 
this purpose, such as fiTQun (56} }5Z] (private communi¬ 
cation, Michael Wilking) and MS-fit [SB] (private com¬ 
munication, Michael Smy), if a pure enough sample of 
hits can be obtained. 

To exploit the correlations between photons in the 
same Cherenkov ring, one must connect the solutions 
from separate PMTs. Consider a vertically downgoing 
muon passing through the center of the detector (the 
considerations generalize). The Cherenkov ring from the 
muon is a circle of uniform intensity moving down the de¬ 
tector walls. Charged particles in a shower are arrayed in 
a small, thin, concave bunch centered on the muon. The 
light from forward electrons adds to that from the muon, 
but the light from each deflected electron makes a tilted 
ellipse that intersects the circle from the muon at only 
two points. When we fit only one PMT hit at a time, it is 
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as if we are azimuthally averaging, turning these ellipses 
into broad circular bands that blend with the light from 
each other and that from the forward particles. 

For showers of small energy, which are the most im¬ 
portant in terms of the frequency of isotope production, 
it should be possible to simultaneously fit the Cherenkov 
rings of all charged particles, or at least the most de¬ 
flected ones. This method may also work for low-energy 
delta rays. 

To exploit the correlations between electrons at the 
same position, one must take into account the angular 
distribution of shower particles. Their light follows hy¬ 
perbolas described by Eq. ^ and clearly visible in Fig. [4] 
We emphasize that the Super-K reconstruction method, 
even with our improvements, does not exploit these hy¬ 
perbolic patterns, which is clearly a missed opportunity. 

For showers of intermediate energy, it should be possi¬ 
ble to fit the portions of the hyperbolas in Fig. [4]that can 
be separated from the light from the muon track. The 
angular distribution of electrons at a given point along 
the muon track determines how the intensity varies along 
the hyperbola. It is probably adequate to focus on the 
integral of this intensity, which reveals the number of 
sufficiently deflected electrons at each position. 

To exploit the correlations between electrons at dif¬ 
ferent positions, one must take into account the longi¬ 
tudinal shower profile. At present, the values of the re¬ 
constructed light from each distance bin along the muon 
track are independent. This allows fluctuations between 
different bins that are larger than the intrinsic ones, 
which likely increases the apparent width of showers. 
There should be a way to enforce consistency between 
the values reconstructed for nearby bins. 

For showers of large energy, it should be possible to 
do forward fitting with a template for a shower of un¬ 
known energy and position along the muon track, assum¬ 
ing something about the average angular distribution of 
shower particles. As shown in Fig. 6 in Ref. 12 a, the in¬ 
trinsic shower fluctuations at 100 GeV are minimal and 
those at 10 GeV are moderate; this method may work to 
even lower energies. 

With these or other new methods based on the physics 
of showers, we are confident that the quality of the re¬ 
constructed light profiles can be significantly improved, 
resolving much smaller showers and multiple showers per 
muon. This will allow much sharper cuts on spallation 
isotopes. 

IV. MUON PROFILE LIKELIHOOD 

Our goal is to improve spallation background rejection 
in Super-K, i.e., the separation of spallation decays from 
neutrino signal events. Currently, Super-K uses spalla¬ 
tion likelihood functions that take variables describing a 
candidate signal event and return a probability that it is 
a spallation decay. 

So far our discussions have been within the framework 


of the Super-K likelihood function for their DSNB anal¬ 
ysis 1, which is based on finding the peak position of a 
muon light profile. Our work in previous sections shows 
how to measure this peak better. 

In this section, we propose a new framework. We build 
a spallation likelihood function based on our faithfully 
reconstructed shower profiles. We quantify its improve¬ 
ment to Super-K spallation cuts. 

There is an important distinction between constructing 
a likelihood function and applying it. When construct¬ 
ing a likelihood, one always knows which primary muon 
made a particular spallation isotope, whereas when ap¬ 
plying a likelihood, one does not know which muon to 
associate a particular event with. We explain the first 
part in Sec. |IV A| and the second in Sec. |IVB| 

A. Spallation Likelihood Functions 

In Super-K, solar neutrinos have a low event rate. In¬ 
trinsic radioactivity backgrounds dominate at low energy 
(< 6 MeV); spallation backgrounds dominate at high en¬ 
ergy (6-18 MeV). Both neutrino events and radioactive 
backgrounds are uncorrelated with cosmic-ray muons, 
and we refer to them as random events. 

A Super-K spallation likelihood C(C , M) evaluates 
how likely it is that a candidate event ( C ) is correlated 
with a muon (M). The larger C, the more likely that this 
C is made by this M, i.e., is a spallation decay. Other¬ 
wise, it is likely an uncorrelated random event. (Below, 
we directly define likelihood functions; Super-K analyses 
use the logarithm of the likelihood, which is equivalent.) 

A good likelihood function reflects the physics of how 
muons make spallation isotopes. (Though well motivated 
on general grounds, the Super-K likelihood functions are 
empirical.) There are several steps to build this func¬ 
tion. First, one picks variables describing a candidate 
event that are statistically different for spallation decays 
and random events. Some obvious choices are the dif¬ 
ferences in time and transverse position between a can¬ 
didate event and its parent muon as well as the muon 
energy loss. A basic assumption that Super-K adopts, 
which we keep, is that these variables are independent, 
i.e., that the likelihood function can be factorized. 

Second, one selects a spallation decay sample along 
with their parent muons. In simulation, this is easy. In 
practice, Super-K selects candidate events that are close 
to muons in time and space, and that have high energy 
(to avoid radioactive backgrounds). This is sufficient to 
select an almost pure spallation sample, due to their high 
rate. 

Third, one selects a random sample with uncorrelated 
muons. In simulation, we simply produce candidate 
events that are uniform in time and space, and randomly 
pair them with muons. In practice, Super-K pairs candi¬ 
date events with muons that follow, instead of precede, 
candidate events in time. 

Finally, one builds every component in a likelihood 
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function. For each variable, one measures the distribu¬ 
tions of this variable for the spallation event sample and 
the random event sample. Then, each likelihood com¬ 
ponent is the ratio of the distributions of the spallation 
sample relative to the random sample. 

The first likelihood function that Super-K developed, 
which is still used for solar neutrino analyses DU, is 

£flat = F t {t) ■ Fi(L trans) F q (Q res )• (ii) 


Here t and A tra ns are the time difference and the trans¬ 
verse distance between a spallation decay and a muon 
track, and Q res is the radiative energy loss of the 
muon (measured by subtracting the light of a minimum- 
ionization muon from the total). We can understand the 
behaviors of these components based on shower physics. 
F t decreases due to the exponential decays of spallation 
isotopes. Fi decreases due to the exponential decrease 
in secondary particle, and thus spallation isotope, den¬ 
sity away from the muon track. F q increases due to ex¬ 
cess energy loss producing more secondary particles and 
spallation isotopes. Their functional forms are given in 
Refs. EM. 


The likelihood function in Eq. (11) does not directly 
reflect shower physics. The Q res variable includes energy 
loss from all showers and low-energy delta rays along a 
muon track. It can be close to the energy of the largest 
shower if that shower is very energetic, but most com¬ 
monly it sums over comparable small showers and low- 
energy delta rays. In addition, this likelihood does not 
include shower position information. 

A new likelihood function Super-K recently developed, 
which is applied to the DSNB analysis 0, is 


£peak = F t (t) ■ Fi(L trans ) ■ F q (Q peak ) ■ F/(L long ). (12) 



FIG. 7. Comparison of the Li ong or 2 component of the three 
different likelihood functions, normalized to equal area. The 
y axis unit is the Cherenkov light intensity, but we also use 
it for likelihoods because they have arbitrary normalization. 
The gray shaded region is the real muon light profile, with 
muon light subtracted. (We cut off the y axis at 5, but the 
light intensity at 23 m goes to 8.) The red star at 8 m indicates 
the position of a spallation event. 


To fully utilize the information about showers in the 
reconstructed muon light profiles, we propose a new like¬ 
lihood, 


Here Qpeak is the total light in the central 4.5 m of the 
reconstructed peak and Ti ong is the longitudinal distance 
(along the muon track) between the peak and the can¬ 
didate event. F[ deceases with respect to the absolute 
value of Li ong because spallation isotopes are most fre¬ 
quently produced in the biggest showers. F q increases be¬ 
cause secondary particle path lengths, and hence spalla¬ 
tion production, are proportional to shower energy. Their 
functional forms are given in Ref. 


-^shower — Ff (f) • Fi(Atrans) ' F z (z ). 


(13) 


The likelihood function in Eq. (12) improved upon that 


in Eq. (11) and reflects some shower physics, though this 


was not recognized as the reason. The variable Q pea k 
would be a good measure of the shower energy if the 
reconstruction were perfect, but it is not accurate, as we 
explained earlier. This likelihood keeps the peak position 
of the biggest shower from reconstruction but discards 
the shape of the shower and smaller showers. 

It is easy to see how our work on improving muon 
light profile reconstruction improves the efficiency of this 
likelihood. First, we can measure the true shower en¬ 
ergy without overcounting or undercounting problems 
(see Sec. Ill C). Second, the functional form of F[ gets 
sharper (see Fig. 12 of Ref. [2HI), which more clearly sep¬ 
arates spallation decays from random events. 


Here 2 is the position along the muon track from where it 
enters the ID. For a muon with a faithfully reconstructed 
shower profile, F z (z ) is the shower profile intensity at 
position z. We emphasize that this shower profile is the 
one shown in Fig. [6] with Improvements 1 + 2 + 3, which 
excludes the light from a minimum-ionizing muon and 
some low-energy delta rays. 


The likelihood function in Eq. (13) fully incorporates 
information about showers. Spallation production at 
each position is roughly proportional to the local sec¬ 
ondary particle path length, which is roughly propor¬ 
tional to the local light intensity. Random events have 
a flat position distribution along the muon track. F z di¬ 
rectly reflects how the probability of spallation produc¬ 
tion varies along a particular muon track. 

Figure 7] illustrates the term containing Ai ong or z in 
these likelihood functions for an example muon. The real 
muon light profile is shown in gray. There are two show¬ 
ers. Most frequently, the isotope would be associated 
with the larger one. However, in this case, the isotope is 
produced in the smaller shower. 
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To emphasize the shape differences of these three like¬ 
lihood functions, we normalize them to the same area. 

at cuts background equally everywhere along the muon 
track. It cuts events in shower regions too weakly and 
nonshower regions too strongly. £ pea k correctly picks out 
the biggest shower along the muon track and cuts events 
in that region with more weight, although its shape does 
not trace this shower perfectly because the likelihood is 
from an average shower profile. Further, when there are 
multiple showers, £ pea k cuts events in the small-shower 
regions too weakly. For £ s hower> the likelihood function 
is the reconstructed shower profile itself. As shown in 
this example, it not only traces the large shower well, 
but it also picks out the small shower. It cuts spallation 
isotopes with weight proportional to the local shower in¬ 
tensity, which is close to optimal. 

B. Efficiency Improvements 

We now quantify the improvement of our shower likeli¬ 
hood function over the Super-K likelihood functions. To 
do so, we need to explain how to apply a cut, i.e., de¬ 
cide whether to discard a candidate event as a spallation 
decay on the basis of a likelihood test. 

When applying a cut to a candidate event, one does 
not know which muon, if any, is correlated with it. This 
is different from building a likelihood, where every spal¬ 
lation decay is paired with its parent muon and every 
random event is paired with an uncorrelated muon. 

The first step of applying a cut is thus to build an event 
likelihood £c(C) that returns the likelihood of an event 
being a spallation decay (from any muon). This is done 
by taking a likelihood function C{C , M) and marginaliz¬ 
ing over all muons {Mi} that are possibly correlated with 
this candidate, which in practice are muons in the previ¬ 
ous 100 s (~ 200 muons). One calculates = C(C,Mi) 
for each muon. The maximum value of Ci is then as¬ 
signed to this candidate event as its event likelihood: 

C c (C) = max.C(C,Mi). (14) 

i 

Second, one obtains a spallation decay event sample 
and a random event sample, and finds the distributions of 
Cc for both samples. The methods to get event samples 
are as described in Sec. |IV A[ except that one discards 
the information about the candidate-muon correlations. 
Then, one calculates the distributions of the event likeli¬ 
hoods Cc- 

For a specific likelihood function, the distributions of 
Cc for the spallation and random samples both have a 
bump, and drop off at small and large values. By design, 
the average Cc value for the spallation sample is larger 
than that for the random sample. However, the two dis¬ 
tributions have significant overlap, which is why it is dif¬ 
ficult to categorize a candidate as a spallation decay or 
a random event. Because of the arbitrary normalization 
of C{C, M), one should not compare Cc distributions for 
different likelihood functions. 



FIG. 8. Efficiency comparison for different likelihood func¬ 
tions. The x axis shows the fraction of spallation background 
events remaining after cuts. The y axis shows the fraction of 
signal events rejected by the cuts. 


Third, one chooses a value, £ C ut, that best separates 
the Cc distributions from the two samples. The choice of 
£ C ut determines its effects on the signal and backgrounds, 
characterized by deadtime and cut efficiency. Deadtime 
is the fraction of random events with Cc > £ C ut , which 
defines the signal loss. Cut efficiency is the fraction of 
spallation decays with Cc > £ C ut, which describes the 
background rejection. Hence, we want to minimize dead¬ 
time while we maximize cut efficiency. For too small a 
value of £ cu t, the deadtime would be unacceptably high. 
For too large a value of £ cut , the cut efficiency would be 
unacceptably low. An optimal value £ cut must be cho¬ 
sen to maximize the signal detection significance. For 
Super-K flat likelihood function, the cut efficiency and 
deadtime for the optimal £ cut are about 90% (10% back¬ 
ground remaining) and 20% for their solar neutrino anal¬ 
yses [131111. 

Finally, one applies the cut to the real data sample, 
rejecting events with Cc > C cut . 

Now we can compare the efficiencies of different like¬ 
lihood functions. We vary the £ cut value for each likeli¬ 
hood, obtaining pairs of deadtime and cut efficiency val¬ 
ues. The optimal £ cut value would correspond to specific 
values. 

To better separate the factors that contribute to the 
differences between likelihood functions, we make some 
simplifications. First, to show the maximum improve¬ 
ment possible due to better reconstruction methods, we 
take the real (simulated) muon light profiles instead of 
the reconstructed ones. Second, to fairly compare the 









14 


difference between the peak and shower likelihood func¬ 
tions, we add an improved peak likelihood that we ex¬ 
plain below. Third, we include only single throughgoing 
muons. Lastly, our spallation samples have only spalla¬ 
tion decays, whereas in Super-K analyses there are some 
random events. Despite these simplifications, our results 
for the flat likelihood are consistent with the Super-K 
measurements and the differences among different likeli¬ 
hoods should be reasonably accurate. 

Figure [8] shows the deadtime and cut efficiency for 
the three likelihood functions. To emphasize the im¬ 
provement, we show background remaining, which is 
(1 — cut efficiency), on the x axis. The flat [Eq. 01- 
peak [Eq. and shower [Eq. fll3| ] likelihoods take 

the functional forms we defined, with components taken 
from Super-K measurements or our definitions. To take 
into account the fact that we use real muon light profiles, 
we show the peak improved likelihood, which is based 
on the Super-K peak likelihood formula, but we adjust 
F[(L\ ong ) (both shown in Fig. 12 in Ref. [20]), assuming 
that the peak position can be measured perfectly. 

To make the comparisons specific, we compare likeli¬ 
hoods at fixed deadtime (~ 20%), which focuses on back¬ 
ground reduction. (One could also compare at fixed cut 
efficiency, which focuses on signal gain.) 

Going from the flat likelihood to the Super-K peak like¬ 
lihood, the background remaining decreases from 0.12 to 
0.09. This shows how the method of Ref. 0 used for 
the DSNB analysis could benefit Super-K solar studies 
by focusing cuts on regions where the muon light pro¬ 
file peaks. With a better peak localization, as our tech¬ 
niques could provide, the improvement would be to 0.07. 
Finally, with our new shower likelihood, the complete im¬ 
provement would be from 0.12 to 0.05. Thus, it should 
be possible to reduce backgrounds by more than a factor 
of 2 with the results we present in this paper. This can 
be combined with other cuts we have suggested or will 
present in forthcoming papers. 

Our proposed new shower likelihood function, based 
on better reconstructed muon light profiles, could sub¬ 
stantially reduce backgrounds. It takes variables directly 
from each muon light profile, so it should be easy to im¬ 
plement. 

In addition to the single throughgoing muons we con¬ 
sider, there are three other classes of muons identified by 
Super-K that are subdominant but relevant ;9]. Stop¬ 
ping muons only make isotopes when fi~ undergo nu¬ 
clear capture [20;, 40], and we discuss cuts in Ref. [20] . 
Corner-clipping muons are just a category of throughgo¬ 
ing muons where reconstruction is more difficult. Iso¬ 
topes will be produced in the FV only when there is 
a shower that enters the ID (excess light in the outer 
detector may help identify large showers), because the 
lateral extent of a shower (< 1 nr) and thus of isotope 
production is less than the thickness of the ID-FV shield¬ 
ing (2 m) [19[ [20] ■ For multiple muons, also known as 
muon bundles, pairs of muons produced in the same at¬ 
mospheric shower are the most common case 


Higher-multiplicity events can be cut aggressively with¬ 
out appreciable deadtime. Our reconstruction method 
could be adapted to deal with pairs, treating them to¬ 
gether when the separation is < 1 m and singly when 
it is larger, along with straightforward adjustments for 
the amount of light and number of showers expected. In 
summary, we see no barriers to adapting our methods to 
implement a complete background-rejection program in 
Super-K. 


V. CONCLUSIONS 

Muon-induced spallation backgrounds are a major 
source of background for low-energy neutrino detection. 
Nevertheless, a complete picture of how these spallation 
isotopes are produced and a strong enough background 
rejection method have been lacking. We are conducting 
a series of studies intended to provide a comprehensive 
theoretical understanding of how muons make spallation 
isotopes and to propose better ways to reject them. 

In our previous papers [19j [20|, we found that almost 
all spallation isotopes are produced by secondary parti¬ 
cles, and that almost all secondary particles are made 
in muon-induced showers. Our calculations agree with 
Super-K measurements on the total spallation yield and 
other data. We also explained an empirical cut that 
Super-K developed for their DSNB analysis [S]. However, 
one discrepancy remained: The Super-K reconstructed 
muon light profiles show prominent bump features, which 
are grossly inconsistent with shower physics. 

In this paper, we show that the observed bump features 
are indeed caused by muon-induced showers. The reason 
that they look too wide and short compared to showers 
can be traced back to the Super-K muon light profile 
reconstruction method. 

We suggest ways to improve the Super-K reconstruc¬ 
tion method. By measuring the position and time for 
every PMT hit, Super-K solves a quadratic equation for 
the emission position of the photon along a muon track. 
However, due to the electron deflection in showers and de¬ 
tector resolution, the quadratic equation could have zero, 
one, or two solutions for one PMT hit, and the solutions 
could be shifted from their true values. We propose ways 
to improve this by picking out one solution when there 
are two, recovering one solution when there seem to be 
zero, and focusing on the PMT hits that give solutions 
closest to the true values. We show that our improve¬ 
ments could lead to almost perfectly reconstructed muon 
light profiles. 

We then propose a new spallation likelihood function, 
based on better reconstructed methods, that fully ex¬ 
ploits the information contained in muon light profiles. 
We demonstrated that it, combined with a better recon¬ 
struction method, can reduce the remaining background 
by a factor of 2 compared to the Super-K DSNB analysis 
cut, and even more compared to their solar analysis cut. 

Our results could be easily adopted by Super-K for 
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their solar neutrino and DSNB analyses. The background 
rejection improvement could be especially dramatic for 
solar neutrinos, where the current cut does not even take 
advantage of the muon light profile peaks, much less the 
full understanding of shower physics. 

The techniques we developed will benefit other neu¬ 
trino experiments. Our results have immediate applica¬ 
tions for other water Cherenkov detectors, e.g., Hyper- 
Kamiokande [44] , which will be shallower than Super-K. 
And, because our reconstruction method does not de¬ 
pend on the geometry of the PMTs, it could be applied to 
muon reconstruction in high-energy neutrino telescopes 
like IceCube |46j . where fluctuations in shower energy 
along muon tracks are used to estimate muon energy. Fi¬ 
nally, our reconstruction technique does not depend on 


the direction of the light, so our results could be adapted 
for scintillator detectors [see Eq. ([8]) for isotropic light 
emission], especially in large-scale next-generation detec¬ 
tors such as JUNO [35] , 
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Appendix: Additional reconstruction examples 

Here we show examples of muon light profile recon¬ 
struction for a variety of other cases. 

Figure [9] shows a muon event with only small showers. 

There are two showers, at ~ 2 and 16 m. The Super-K 


profile using the Super-K method does not reveal either 
shower, and is mostly noise. Our reconstructed profile, 
however, successfully picks out both showers. 

Figure [T0| shows a muon event with moderate but com¬ 
parable showers, at ~ 12 and 23 m (and smaller ones at 
~ 2 m and 30 m). The Super-K reconstructed profile is 
reasonable for one shower, poor for the others, and has 
a false shower in the middle. Our reconstructed profile 
reconstructs all four showers clearly. 

Figure [lT] shows a muon event with a hadronic shower. 
This is a relatively clean event. Both profiles have the 
correct peak position, and trace the shower shape. The 
Super-K profile is more smeared out. Our reconstructed 
shower is not as sharp as large in other examples due 
to the larger lateral displacement of charged particles in 
hadronic showers. 

In summary, these examples demonstrate the better 
performance of our reconstruction for small, multiple and 
hadronic showers. 
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FIG. 9. Shower reconstruction example—only small showers. Left panel: The gray shaded shape is the real (simulated) 
light profile. The black line is the result of the reconstruction using the Super-K method. The blue line is the result of our 
improvements. Right panel: The time and position pattern of Cherenkov light on the ID walls for this muon event. 
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FIG. 10. Shower reconstruction example—two comparable showers. Descriptions as in Fig. [9] The only difference is that we 
use a bin size of 1 m in the left panel to reduce fluctuations. 
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FIG. 11. Shower reconstruction example—hadronic showers. Descriptions as in Fig. [9] 











































































